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ABSTRACT 

We explore the possibility of detecting radio emission in the cosmic web by analyz- 
ing shock waves in the MareNostrum cosmological simulation. This requires a careful 
calibration of shock finding algorithms in Smoothcd-Particle Hydrodynamics simula- 
tions, which we present here. Moreover, we identify the elements of the cosmic web, 
namely voids, walls, filaments and clusters with the use of the Spine Web technique, a 
procedure that classifies the structure in terms of its topology. Thus, we are able to 
study the Mach number distribution as a function of its environment. We find that the 
median Mach number, for clusters is M. dusters ~ 1-8, for filaments is .M n i am0 nts ~ 6.2, 
for walls is -MwaUs ~ 7.5, and for voids is Avoids ~ 18. We then estimate the radio 
emission in the cosmic web using the formalism derived in Hoeft & Briiggen (2007). 
We also find that in order to match our simulations with observational data (e.g., 
NVSS radio relic luminosity function), a fraction of energy dissipated at the shock of 
£ = 0.0005 is needed, in contrast with the £ c = 0.005 proposed by Hoeft et al. (2008). 
We find that 41% of clusters with M > 1O 14 M host diffuse radio emission in the 
form of radio relics. Moreover, we predict that the radio flux from filaments should be 
Si50MHz ~ 0.12 /xJy at a frequency of 150 MHz. 

Key words: cosmology: theory - large-scale structure of the Universe - hydrody- 
namics - methods: numerical - radiation mechanisms: nonthermal - shock waves 



> ' 1 INTRODUCTION 



In the standard theory of formation history, matter evolved 
from small perturbations in the primordial density field into 
a complex structure of sheets and filaments with galaxy clus- 
ters at the intersections of this filamentary structure. Galaxy 
surveys, such as the 2dF-GRS ( 2 degree Field Galaxy Red- 
shift Survey, IColless et al]|2003) and the SDSS (Sloan Digi- 
tal Sky Survey. iTegmark et alj 20041 ) . have revealed a com- 
plex network of filamentary natu re, which has bec ome nec- 
essary known as the cosmic web (Bond et al. 19961). 

In the evolution of the cosmic web, baryons in the dif- 
fuse intergalactic medium accelerate towards dark matter 
halos under the growing influence of gravity and go through 
shocks that heat them. These cosmological shocks are a 
ubiquitous consequence of cosmic structure formation. They 
are tracers of the large-scale structure and contain infor- 
mation about its dynamical history. Gas in walls and fila- 
ments follows the gravitational potential towards clusters of 
galaxies, colliding with the intracluster medium (ICM) at 
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speeds of > 1000 km s _1 . While cosmological shocks heat 
the ICM to temperatures of ~ 10 7 — 10 8 K, filaments are 
heated to temperatures of 10 5 — 10 7 K, making their detec- 
tion at present time challenging, since at this temperature 
there are hardly any emission lines. Therefore, the study of 
shocks may provide an independent and complementary way 
to study the Warm-Hot Intergalactic Medium (WHIM), the 
low density intergalactic medium of the cosmic web that is 
believed to host the majority of the fraction of the baryon 
density. 

Some cosmological shocks are associated with diffuse 
radio emission caused by electrons a ccelerated to rela- 
tivistic velocities by Fermi acceleration (|Enfilin et alj 1 19981 : 
Rocttigcr et al . 1999: lMiniati et alj|200lh . This diffuse radio 
emission, without galaxy counterpart, is usually divided into 
two classes, namely radio halos and radio relics. Radio ha- 
los are unpolarized and have diffuse morphologies that are 
similar to those of the thermal X-ray emission of the clus- 
ter gas (|Giovannini et al1l2006l ). They are usually found in 
the center of cluster s with significant substructure (see e.g., 
ICassano et"afll201Ch . 

Radio relics, on the other hand, are typically located 
near the periphery of the cluster. They have been observed 
in several merging galaxy clusters. They often exhibit sharp 
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emission edges and many of them show strong radio polar- 
ization. Once accelerated, the electrons are short-lived be- 
cause of the inverse Compton scattering and synchrotron 
energy losse s, and their spectrum rapidly steepens from the 
shock edge l|Giacintucci et al1l2008l ; Ivan Weeren et aLll2009l ) 
The synchrotron nature of this radio emission indicates the 
pres ence of cluster magnetic fie lds of the order of ~ 0.1 — 1 
MG jFeretti fc Giovanninill2008h . 

Cosmic shocks in large-scale structure s have been in- 
vesti ga ted in a number of semi-ana ly tical jGabici fc Blasil 



20031; Berrington 



Dermerl 120031 ; iKeshet et al.1 120031; 
Meli fc Biermannll2006l) and numerical work (fMiniati etakl 
200dl200ll;lMiniatill2q02];lRvu et al.ll2003l; iKang et all 
Pfrommer et al.l l200fl 120071. 12003: 



Jubclgas et al. 



Hocft et al. 2008; Skillman ct al. 2008;: Battaglia et al. 



2005 



2008 



2009; 



Vazza et al.l 120091 ; Iskillman et al.l 120101 ) . Some studies have 



focused on the non-thermal emission from galaxy clusters 
by modeling discretized cosmic ray (CR) energy spectra. 
A series of papers explored the dynamical impact of cos- 
mic ray (CR) pr otons on hydrodynamics in cosmological 



SPH simulations (IPfrommer et al-feoOrj ; lEnfilin et al]|2007l ; 
ljubelgas et al.l I2008T ). O bservat ionally, detecting shocks 
waves in large-scale structures is still challenging, since 
they usually develop in the external regions of galaxy clus- 
ters, where the X-ray emission is faint. However, a few 
merger shocks h ave been detected in n e arby X-ray bright 
galaxy clusters dMarkevitch et~ai1 120051 ; iMarkevitchl 120061; 
ISolovveva et al.ll2008l ) and may be possible associated with a 
single or doub le radio relics discovered in a number of galax 



2008). 



clusters (e.g., RSttgcring ct al. 1997; Markcvitch ct al 
Bagc hi et alj 120061 ; iBonafede et alj 120091 ; iGiacintucci et al 



alaxy 
2005; 



In numerical studies, iMiniati et al.l (|2000r i studied the 
quantitative properties of large-scale shocks produced by gas 
during the formation of cosmic structures by means of hy- 
drodynamical simulations. They showed that shocks form 
abundantly in the course of structure formation and their 
topology is very complex and highly connected. They also 
stated that considering the large size and long lifetime of 
shocks, they are potent i ally in teresting s ites for cosmic-ray 
acceleration. iRvu et all l|2003h (see also IKang et"al1 l2005l l 
found that shocks form around sheets, filaments and knots 
of mass distribution when the gas in void regions accretes 
onto them. 

|pfr ommer et al.l l|2006l ) studied the properties of struc- 
ture formation shock waves in cosmological SPH simula- 
tions, allowing them to study their role in the thermalization 
of the plasma as well as for the acceleration of relativistic 
CRs through diffusive shock acceleration. They find that 
most of the energy is dissipated in weak internal shocks, 
with Mach number M ~ 2. On the other hand, collapsed 
cosmological structures are surrounded by external shocks 
with much higher Mach numbers (Ai ~ 100), but they play 
only a minor r o le in the energy balance of thermalization. 
Skill man et al.l |2008l 'l computed the production of CRs in 
a cosmological simulation volume, finding that CRs are dy- 
namically important in galaxy clusters. They also found that 
shocks with low Mach number typically trace mergers and 
complex flows, while those with mild and high Mach num- 
ber (A4 > 5) g enerally follow acc retion onto filaments. In 
a similar study, IVazza et al.l (|2009h analyzed the properties 
of large-scale shocks. They find that the bulk of the energy 



in galaxy clusters is dissipated at weak shocks, with Mach 
numbers A4 ~ 1.5, although slightly stronger shocks are 
found in the external regions of merging clusters. 

IPfrommer et all i^POST ) used GADGET simulations of a 
sample of galaxy clusters, implementing a formalism for CR 
physics on top of radiative hydrodynamics. They modeled 
relativistic electrons that are accelerated at cosmic forma- 
tion shocks and produced in hadronic interactions of CRs 
with protons of the ICM. They found that the radio emis- 
sion in clusters is dominated by secondary electrons. Only 
at the location of strong shock s the contribution of pr imary 
electrons may dominate. Later. [Battaglia et al.l (|2009h stud- 
ied the radio emission for a sample of 10 galaxy clusters 
extracted from zoomed cosmological simulations. They de- 
termined the radio luminosity, spectral index and Faraday 
rotation measure for regions of the relics and concluded that 
upcoming radio telescopes, such as LOFAR, MWA, LWA 
and SKA will discover a substantially larger sample of ra- 
dio relics. SKA should probe the macroscopic parameters of 
plasma physics in clusters. 

iHoeft et alj (|2008l ) investigated the diffuse radio emis- 
sion from clusters in cosmological SPH simulations. They 
found that the maximum diffuse radio emission in clusters 
depends strongly on their X-ray temperature. They also 
found that the so-called accretion shocks cause only very 
little radio emission. They conclude that a moderate effi- 
ciency of shock acceleration, namely £ c = 0.005, together 
with moderate magnetic field, namely 0.07-0.8 /j,G, in the 
region of radio relic are sufficient to repro duce the num- 
ber d ensity and luminosity of radio relics. ISkillman et al.1 
(2010), using AMR simulations, also studied radio emission 
in galaxy clusters. They investigated scaling relations be- 
tween cluster parameters such as synchrotron power, mass 
and X-ray luminosity. 

In this paper we follow the approach of IHoeft et al.1 
(2008) but, instead of focusing on galaxy clusters, we inves- 
tigate whether it is possible to detect radio emission in the 
entire cosmic web. Due to the low temperature of the accret- 
ing gas, the Mach number of external sho cks - filaments - is 
high, extending up to M ~ 100 or higher (|Rvu et al.l 1 20031 ). 
Therefore, the total radio power generated by infall along 
filaments will be significantly weaker than infall into a clus- 
ter. We will try to determine this total radio power. To this 
end, we use the approach for detecting and identifying shock 
fronts and their c haracteristic Mach number developed by 
Hocft et al.l (|2008r ). This shock finder is applied to a large 
iV-body/SPH simulation. For com puting the radio emission , 
we follow the model elaborated in IHoeft fc Briiggenl |2007). 
They assumed that electrons are accelerated by diffuse shock 
acceleration and cool subsequently by synchrotron and in- 
verse Compton losses. Consequently, the radio emission can 
be expressed as a function of downstream plasma properties, 
Mach number and surface area of the shock front. In order 
to correctly assign a surface area to the SPH particles, we 
calibrate the radio emission using shock tubes experiments. 
Of primary importance is the proper identification of the 
different structures that reside in the cosmic web. To this 
end, w e apply the Spine Web technique (|Arag6n-Calvo et al.l 
l2010al . see also lAragon-Calvo et al.ll2010d ) to the entire sim- 
ulation. The Spine Web procedure correctly identifies voids, 
walls, filaments and clusters present in the cosmic web. It 
is a powerful method that deals with the topology of the 
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underlying density field of the cosmic web. This allows us 
to study shock properties and radio emission in the differ- 
ent structures that constitutes the cosmic web. This goes 
beyond the traditionall y internal v e rsus e xternal shock clas- 
sification suggested by iRvu et al.l |2003l ) and the method 
of characterizing shocks b y their pre-shock over density and 
temperature proposed by ISkillman et all ([2008) . Applying 
the radio emission model to the simulation leads to radio- 
loud shock fronts which we are able to locate within the 
cosmic web. 

This paper is organized as follows. In Section[2]we sum- 
marize the characteristics of the A-body/SPH simulation. 
Section [3] describes the Spineweb method. In Section [4] we 
give a brief description of the physics involving non-radiative 
shocks and depict how we find shocks in SPH simulations. 
In order to calibrate the radio emission model, we use shock 
tube tests which are described in Section [6] The radio emis- 
sion model is sensitive to the shock surface area. Using the 
results of the shock tubes, in Section 16.11 we described the 
method used in order to precisely determine the area of the 
shock front. In Section[7]we present our results for the shock 
fronts and radio emission in the MareNostrum Universe. In 
Section [8] we summarize our findings. 



2 THE COMPUTER SIMULATION 

In order to study cosmological shock waves in the large-scale 
structure o f the Universe, we use the MareNostrum Universe 
simulation ijGottlober fc Yepesll2007l l , one of the largest cos- 
mological simulations available. It assumes a standard flat 
ACDM Universe with cosmological parameters f2 m ,o = 0.3, 
Q,b,o = 0.045, fU,o = 0.7 and h = 0.7, where the Hub- 
ble parameter is given by Ho = 100ft km s _1 Mpc _1 . The 
normalization of the power spectrum is erg = 0.9. The sim- 
ulation box has a side length of 500/i _1 Mpc and contains 
1024 3 gas and 1024 3 dark matter particles, with masses of 
m gaa = 1.45 x 10 9 /i _1 M© and m dm = 8.24 x 10 9 /i _1 M Q . 
The initial conditions are followed from a redshift of z — 40 
until the present time (z = 0) using the massively parallel 
tree TV-Body /smoothed par ticle hydrodynamics (SPH) code 
GADGET-2 ilSpringellbOol ). The results presented here are 
at present time (z = 0). Radiative processes or star for- 
mation are not included in the simulation. The Plummer- 
equivalent softening was set at epi — 15 /i _1 kpc in comov- 
ing units, and the SPH smoothing length was set to the 
distance to the 40th nearest neighbor of each SPH parti- 
cle. Smoothing scales are not allowed to be smaller than the 
gravitational softening of the gas particles. 

To extract the groups and galaxy clusters pre sent in 
the simulation, we use HOP l|Eisenstein fc Hutll 19981 ). HOP 
assigns a density to every particle by smoothing the density 
field with a spline cubic kernel using the n nearest neighbors 
of a given particle. Particles are then linked by associating 
each particle to the densest particle from the list of its closest 
neighbor. This process is repeated until it reaches a particle 
that is its own highest density neighbor. All particles linked 
to a local density maximum are identified as a group. At 
this stage, no distinction between a high-density region and 
its surrounding has been make. To identify halos above a 
certain density threshold, a regrouping merging procedure 
is performed. The code first includes only particles that are 



above some density threshold. Then, it merges all groups for 
which the boundary density between them exceeds a certain 
density value. Finally, all groups identified must have one 
particle that exceeds a density peak to be accepted as an 
independent group. We associate this density peak with the 
virial density needed for a spherical region to be in virial 
equilibrium, A V ir- The value of A v i r is obtained from the 
solution to the collapse of a spherical top-hat perturbation 
under the assumption that the object has just virialized. 
For the cosmology described here, at present time its value 
is A v i r ~ 337. 

In using HOP, we use the dark matter particles of the 
MareNostrum simulation. Once we identified the virialized 
halos of the sample, we take the position of the densest dark 
matter particle as a first estimate of the center of mass. Sub- 
sequently, we add the gas particles, grow a sphere around 
the center of mass and begin iterating, shrinking the sphere 
around the (new) center of mass until we reach a minimum 
of 50 particles. This ensures a correct identification of the 
center of mass taking into account gas and dark matter par- 
ticles. Once we have the final center of mass, we grow a 
sphere around it that encloses a value of A v i r ~ 337. Given 
that we want to study radio emission in the cosmic web, 
the nodes of the filaments will be consider as galaxy clusters 
with masses M > 1O 14 /i _1 M . We find that the MareNos- 
trum Universe contains 3865 clusters of galaxies, the most 
massive with a mass of M = 2.5 x 10 15 /i _1 Mq and virial 
radius R V1I = 2.75/i _1 Mpc. 



3 THE SPINEWEB PROCEDURE 

The characterization of the Large-Scale Structure was done 
usi ng the SpineWeb method (for a t horough explanation 



_ Aragon-Calvo et al.ll2~010al . see also lAragon-Calvo et alj 
l2010ch . This method classifies the Cosmic Web into four ba- 
sic morphological and dynamical components: voids, walls, 
filaments and clusters. T he SpineWeb method is base d on 
the Watershed transform l|Beucher fc Lantueioul|[l979h and 
its cosm ological implement ation: the Watershed Void Finder 
(WVF) dPlaten et alj|2007h . The idea behind the WVF is to 
segment the density field into individual basins by "flood- 
ing" it (in analogy to a landscape). The SpineWeb method 
takes this analogy and goes one step further by directly re- 
lating the topology of the density field and the boundaries 
between basins by observing the simple relation between 
number of adjacent basins and topology: walls correspond 
to the boundary between two voids. Filaments+clusters are 
found at the intersection of three or more voids, which trans- 
lates to the intersection of three or more walls: 



= 1, void 
= 2, wall 
> 3, spine 
(filament + clusters) 



(1) 



The resulting morphological characterization is 
parameter- free. The analysis presented here is also based on 
the most recent implementation of the SpineWeb method 
which extends the original formalism on a multiscale 
hierarchical fashion making it also scale-independent. 

In practice we start our analysis by generating a low res- 
olution version of the original simulation using the averaging 
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4 HYDRO DYNAMICAL SHOCKS AND 
SHOCK FINDER 

In the formation of structures in the Universe, the bulk ve- 
locity of the gas flow often exceeds the local sound speed. 
As a result, the shock surface separates two regions: the 
upstream regime or pre-shock regime, and the downstream 
regime or post-shock regime. The plasma in the upstream 
regime moves with a velocity v u towards the shock front, 
while the downstream plasma departs with velocity ivj. Most 
of the incident kinetic energy flux is converted into thermal 
energy. As the plasma passes through the shock front, mass, 
momentum and energy fluxes are conserved, which is ex- 
pressed in the Rankine-Hugoniot relations: 




Walls and Filaments 



Figure 1. An example of how the SpineWeb procedure divides 
the cosmic web. Shown are surfaces enclosing voxels which are 
identified as belonging to walls (blue, top left frame) and fila- 
ments (red, top right frame) within a cubic region of 50h~ iJVIpc. 
The bottom frame shows the same region with both components 
connected and intertwined. 



procedure described in lKlypin et al.l (|200ll ). The lower par- 
ticle resolution was chosen to correspond to a inter-particle 
separation of approx. 3/i -1 Mpc at the initial conditions. 
With this "linear-regime" filtering we avoid small-scale sub- 
structure at the present time while keeping the characteris- 
tic large-scale anisotropy of the cosmic web. From the par- 
ticle distribution we compute the density field u s ing the 
DTF E method dSchaap fc van de WevgaertJ l200d : ISchaad 
120071 : | van de Weveaert fc Schaad 12009). The field is sam- 
pled at 27 different points inside each voxel and its mean 
is used to estimate the density for that particular voxel. 
From the density field we compute the watershed transform 
and identify the boundaries of the basins (i.e. the water- 
shed transform itself). Finally for each voxel in the water- 
shed transform we apply the criteria outlined in Eqn.Q] The 
complete procedure takes a few minutes on a regular work- 
station. 

The SpineWeb method provides a complete framework 
for the characterization of the cosmic web into its basic mor- 
phological constituents: voids, walls, filaments and clusters. 
It has no free parameters and being fundamentally a topo- 
logical measure it is highly robust against noise in the den- 
sity field. Fig. [T] shows an example of how the SpineWeb 
correctly identifies the different elements of the cosmic web. 
Shown is a full 3D network of filaments (red) and walls 
(blue). It is clearly visible the three-dimensional nature of 
the filament-wall network. Filaments define an intercon- 
nected web and walls fill the spaces in between the filaments 
forming a closed "watertight" network of voids. 



Pd + Pd«d , 
Y +u d + — 



(2) 



puV u = Pd^d • 
Pu + PuV u 
Vl , , Pu 

2. p u l pd 

where p denotes the mass density, P is the pressure and u is 
the specific internal energy, and the velocities are measured 
in the rest-frame of the shock surface. The strength of a 
shock is given by the Mach numbers M, 



(3) 



where c s denotes the upstream sound speed which depends 
on the the specific internal energy by c u = \J 7(7 — l)w u . 
Combining Eqns. [2] and |3j and assuming that the plasma 
obeys the polytropic relation, we can write for the Mach 
number 



r qr 

7 r 



(4) 



where r and q are the compression and entropy ratios. The 
difference between downstream and upstream velocity can 
be written as 



I'd 



Wu = M - 



(5) 



The Mach number estimator, which will be explained in 
the next section, computes the entropy ratio and the ratio 
(v c \ — v u )/ud for each SPH particle in the simulation. In 
order to obtain the Mach number, we tabulate the relations 
between q, r, Ai and J\A{r — l)/r, allowing us to simple read 
the Mach number from a table. 



4.1 Shock finder and Mach number estimate 

In order to calculate the radio emission present in the 
MareNostrum simulation, the first step is to identify shock 
fronts, and from this, derive the Mach number, downstream 
(postshock) temperature and electron density. The method 
we em ploy is the one proposed and used by Irloeft et al.l 
(2008). For a detailed description of the method, we refer 
to the mentioned paper. Here, we will give the basic details. 

The first step consists of computing the entropy gradi- 
ent, VS, for each SPH particle. The entropy gradient gives 
the direction of the shock normal pointing into the down- 
stream direction. We define an associated upstream and 
downstream position, 



xh = Xi + fhhml , 



(6) 
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where I = u,d denotes upstream or downstream, with the 
downstream position in the opposite direction, is the posi- 
tion of the SPH particle i, hi is the smoothing length of par- 
ticle i and nj denotes the shock normal, nj = — VS/|VS|. 

The upstream and downstream velocity are computed, 
together with the internal energy and density, using the 
usual SPH scheme. The upstream velocity is given by 
v U i + v s h = v(x u i) ■ n}, an d the same for the downstream 
velocity, but using the downstream position. v s h is the ve- 
locity of the shock front in the rest-frame of the simulation. 
It is important to notice that the velocity field can also show 
perpendicular components. We compute these components 
in a similar way to the upstream and downstream velocities, 
u* = v(xi ± fhtiiD-i) ■ nf, where the three vectors n 1 , n 2 
and n 3 form an orthonormal base. 

Once this is done, for those particles that belong to a 
shock front we demand that the velocity has to be divergent 
in the direction of the shock normal, i.e., (vd —v n ) > 0. Also, 
the velocity difference in the directions perpendicular to the 
shock normal have to be smaller than that parallel to the 
shock, i.e., |v+ - v~\ < (v d - v u )/2. 

Finally, to compute the Mach number, there are two 
options: it is possible to use the entropy ratio q — Sd§u (S 
is the entropy) or use the ratio (w<j — v u )/c u , where c u is the 
sound speed in the upstream region (see Section^]). In order 
to have a conservative estimate for the Mach number, the 
smaller value of these two options is used. 



5 RADIO EMISSION MODEL 

In this section, we will briefly describe the radio emis- 
sion model. For a co mpl ete descript i on, w e refer to 
iHoeft fc Briiggenl (|2007l ) and lHoeft et~ai1 (|2008l ). 

The first assumption of the model is that radio emis- 
sion is produced by electrons accelerated to ultra-relativistic 
speeds required for synchrotron and inverse Compton radi- 
ation. These electrons are accelerated to a power-law dis- 
tribution that is related to the Mach number from diffusive 
shock acceleration theory. As for the magnetic field in the 
downstream region, it is assumed that, on average, simply 
obeys flux conservation to follow: 



_Bd_ = / n d \2/3 
B Icf ~ Vl0- 4 cm" 3 / 

where B TC f is 0.1/iG and is the number density in the 
downstream region. This is a simple assumption since a de- 
tailed model for the generation and evolution of magnetic 
field s is complex and beyo nd the scope of this paper (see, 
e.g., iBonafede etTafl l2010l for a recent study on magnetic 
fields in the Coma cluster). 

There are several indications that the magnetic eld in- 
tensity decreases going from the center to the periphery of a 
cluster. The exact power with which it varies with density is 
unclear, however. Some simulations suggest a higher power 
than that given by flux conservation which i s attributed 
to the wor k of a fluctuation dyn amo (see e.g. iDolag et al.l 
(2008) and Briigge n et all (|2005 | )). Faraday ro t ation mea- 
surements in the Coma cluster IBonafede et al.l |20ld) give 
a best fit for B oc n 0,5 . Such a power would decrease the 
luminosity of cluster shocks by about 20%. Outside of clus- 
ters there are even fewer observational constraints. Recent 



constraints from TeV 7-ray sources have placed lower limits 
on extragalactic magnetic fields (EGMFs) . The charged par- 
ticles from pair cascades are deflected by EGMFs, thereby 
reducing the observed point-like flux. Dolag et al. (2011) 
have calculated the fluence of 1ES 0229+200 as seen by 
Fermi-LAT for different EGMF. The non- observation of 1ES 
0229+200 by Fermi-LAT suggests that the EGMF fills at 
least 60% of space with fields stronger than 10 -15 G assum- 
ing that the source is stable for at least 10 4 yr. The fields 
outside clusters, however, will have a negligible effect on the 
radio luminosity function since the by far the largest part 
of the radio luminosity originates from shocks in or near 
clusters. 

The radio emissivity, P(e, v bs), which depends on three 
quantities: i) the energy of the electron energy (e), ii) the ob- 
serving frequency (z/ bs), an d iii) the magnetic field strength 
(B), is then given by: 



dP(t/ obs ) 



6.4 x lO^ergs^Hz -1 



Mpc 2 10- 



n e _3 

r cm 



X 0.05 V 1.4GHz/ 

{B/„G) 1+S/2 
(Bcmb/mG) 2 + {B/fiG) 2 



3/2 



7keV, 
9(M,T d ). 



(8) 



where n e is the electron density, £ e is the fraction of energy 
dissipated by the electron at the shock front and s is the 
compression ratio at the shock front. The radiation spectral 
index is related to the electron spectral index by a = (s — 
l)/2. Bqmb has a value of 3.24 = ^iG(l + z) 2 . This magnetic 
field is defined as the magnetic field corresponding to the 
energy density of the CMB. 

As stated before, in order to calculate the radio emis- 
sion, we need the Mach number, M, the downstream tem- 
perature, Td, the electron density, n e , the magnetic field, B, 
and the shock surface area, A. The shock finder described 
in Sect. 14.11 besides locating shock discontinuities in a SPH 
simulation, also provides estimates for the Mach number and 
for the shock normal. This normal vector gives a position 
which is sufficiently downstream to determine Td and n e . 
As for the shock surface area A, each particle represents an 
area of the shock, since a shock front is comprised of particle 
in a SPH simulation. The shock discontinuity is smoothed 
by the size of the SPH kernel, /isph, which contains -/Vsph 
particles, with a corresponding volume of /i|ph- Therefore, 
one particle belonging to the shock front represents a shock 
area of 



A = / A -^ PH 



(9) 



wh ere /a is a normaliz ation constant which was set to be 6.5 
by IHoeft et "all (|2008l ) using shock tube simulations. Here, 
we will correct this factor also using shock tubes, as it was 
found that /a depends strongly on Ai for low Mach number 
shocks (see next section). 



6 CALIBRATING THE RADIO EMISSION 

An important quantity in Eqn. [S] is the shock s urface area, 
which depends on the normalization factor / A . IHoeft et al.l 
(2008) determined this factor using shock tube tests, finding 
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Figure 2. Simulations of the Sod shock tube problem with Mach numbers M = 3, 10 and 1000. The numerical result of the volume- 
averaged hydrodynamical quantities within bins with a spacing equal to the inter-particle separation of the high-density medium is shown 
in open circles while the analytical result is shown as continuous lines. The upper row shows the density p of the gas, the middle row 
shows the velocity v x , while the bottom row depicts the pressure of the gas. 
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0.046 


1000 


1.07 


803.53 


0.050 


1300 


0.63 
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0.051 
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Table 1. Shock tube tests parameters. The columns indicate the 
strength of the shock, i.e., the Mach number M, the pressure in 
the high density region, the fitted Mach number and the standard 
deviation, <r, of the fitted Gaussian. 



that it depends on the Mach number of the shock. The origin 
of this is the following: the Mach number estimator results in 
a distribution of Mach numbers for the particles in the SPH- 
broadened shock discontinuity, i.e., in the periphery of the 
shock front the Mach number is underestimated. For small 
Mach numbers, M < 5 the radio emission depends strongly 
on M, while for large Mach numbers it does not. Therefore, 



a constant /a underestimates the radio emission of shocks 
with M. < 5. Next, we determine the function /a(M) by 
calibrating it with shock tube tests. 

We run several shock tube tests in order to model a 
large range of Mach numbers and, therefore, obtaining a 
prope r area fun ction. We construct 16 standard shock tube 
tests (ISodll 19781) using a 3D- version of the GADGET2 code 
l|Springelll2005l ). We consider an ideal gas with 7 = 5/3, ini- 
tially at rest. The left-half space (x < 2000) is filled with gas 
at unit density, p2 — 1 and pressure P2 = 5 x 10 6 , while the 
right half, x > 2000, is filled with low-density gas, p\ = 0.25 
and variable low-pressure. The value of this low-pressure gas 
has been chosen such that the resulting solutions yield the 
Mach numbers in the range of M = 1.4 — 1800. We set 
up the initial conditions in three dimensions using an irreg- 
ular glass-like distribution. In order to test the differences 
in resolution, we ran the same tests in different resolutions: 
from 12500 particles to 800000 particles. We did not find 
any significant difference between these resolutions. Hence, 
the results presented here are done with shock tubes with 
100000 particles. The particles are contained in a periodic 
box which is longer in the x-direction than in the other two 
dimensions, y and z. The parameters of the different shock 
tube tests are shown in Table [1] 

Fig. [2]shows the profile of gas density, velocity and pres- 
sure in a shock tube calculation where the gas particles ex- 
perience a shock of Mach number M = 10. The simula- 
tion results are presented by circles and the continuous lines 
give the analytic solution. The shock tube simulations were 
run with the same SPH parameters (e.g., artificial viscosity) 
as the MareNostrum Universe in order to compare results. 
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Overall, there is a good agreement between the analytical 
and the numerical solution, with discontinuities resolved in 
about two to three SPH smoothing lengths. There is a char- 
acteristic pressure blip at the contact discontinuity. The ir- 
regularities observed in the density profile may be due to 
the small number of particles in the SPH kernel, namely 
TVsph = 40. 

We calculate the strength of the shock with the Mach 
number estimator. Although the shock tube was run with 
iVspH = 40, the shock finder was set up with -/Vsph = 80 
in order to smooth the fields. To determine the function 
fx{M), we calculate the differential shock surface area in 
logarithmic Mach number bins, dS/d\og M, for each of the 
shock tubes. Fig. [3] shows such distribution for the shock 
tube calculation with M = 10. The surface area was calcu- 
late using Eqn. [9] setting / a = 1 in order to obtain f a (M). 
Visual inspection allow us to draw two observations. Firstly, 
there is a small shift in the value of the theoretical Mach 
number and the one estimated by the shock finder mean- 
ing the shock finder is underestimating the Mach number. 
Secondly, there is a "tail" of shock surface areas to the left 
of the Gaussian-like distribution. This has the following ori- 
gin: the shock finder identifies low-Mach number shocks in 
the periphery of the transition region between the upstream 
and the downstream zone (see Fig. [2} given the nature of 
SPH, i.e., quantities are calculated over a region determined 
by -/Vsph. Hence, the tail has surface area contribution from 
the upstream and downstream region. 

We calculate the shock surface area of all shock tube re- 
alizations. Every shock tube exhibits a similar dS/d log M 
distribution as the one shown in Fig. [3] We fit the distribu- 
tion with a Gaussian centered on the mean Mach number M 
(which should be the theoretical Mach number) and a Heav- 
iside function H multiplied by a constant factor /C which is 
determined as the mean value of the surfaces areas located 
approximately 3a to the left of the Gaussian mean value, 



/(log-M) 



K. ■ H(X) for log log M - 3a 
e -x 2 /2<r 2 for w x > log - 3a 



(10) 



where X — log M. — log M . 

Each shock tube gives us 3 different fitting values: the 
constant /C, the mean of the Gaussian function M and the 
standard deviation of the Gaussian, a. Fig. [4]shows the best 
fits to the mentioned variables. Although we were able to fit 
the constant K, in two segments of Mach numbers, we find it 
is, on average, ~ 10% of the peak value of the Gaussians. For 
M and a, we see there is a good correspondence with the 
theoretical Mach number. We can, therefore, do a linear fit, 
obtaining functions that will deliver a correct Mach number 
M that will be applied to the radio emission model. 



6.1 Calculating the Shock Area 

To correctly calculate the radio emission, we need to correct 
the normalization constant }a{M) (see Eqn. [9}. To this end, 
we find the area within 2a around the mean Mach number 
of each shock tube, S e *_ P - This area should be similar to the 
cross section of the tube, Sthco- We find a relation between 
the Mach number and the ratio Stheo/ S exp : 




10 100 
Mach 



1000 



Ja(M) = 0.541ogX + 5.71. 



(11) 



Figure 3. Differential shock surface area as a function of logarith- 
mic Mach number bins for a shock tube calculation with A4 = 10. 
The differential shock surface area is normalized by the binsizc 
and the cross section of the tube. Solid line represents the bestr 
fit of the Gaussian-like distribution and the constant surface area 
to the left of the Gaussian. 



The fitting obtained in the previous section together 
with Eqn. 1111 allows us to correctly estimate the radio emis- 
sion (Eqn. [5} in the simulation. 

Another important quantity we want to determine from 
the MareNostrum simulation is the shock surface area in a 
given logarithmic Mach number interval. However, as seen 
in the previous section, our shock finder algorithm under- 
estimates high Mach numbers and, due to the SPH nature, 
for high Mach numbers it identifies low shock Mach num- 
bers to particles involved in the shock front. In order to 
correct this, we will use the results of the previous section 
and deconvolve them with the measured dS/dlog M distri- 
bution, i.e., we will deconvolve the shock surface area distri- 
bution calculated from the MareNostrum simulation with a 
deconvolution kernel obtained by constructing the function 
of Eqn. 1101 We do this by using a deconvolution technique 
which is explained in Appendix IA1 

Fig. [5] shows the measured dSj log M distribution ob- 
tained from the MareNostrum simulation (solid line) and the 
"real" distribution obtained with the deconvolution tech- 
nique (filled circles linked by a dotted line), where we have 
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Figure 4. Best fits to the different variables of the Gaussians. 
Upper panel: linear fit of the measured Mach number. Central 
panel: linear fit of the <r of the Gaussian distribution. Bottom 
panel: linear fit to the constant K. left to the Gaussians fits. 



also used the Mach number shift obtained from fitting M 
from the shock tubes. As expected from the shock tube data, 
the real distribution has lower Mach numbers, expressed in 
a smaller shock surface area, while it has more high Mach 
numbers for M > 10. 

The Mach number determination is constructed in a 
way that it results in a conservative estimate even in more 
complex situations. We would like to demonstrate this by a 
colliding shock tube simulation. We have set up a shock tube 
with three initial zones, with densities and internal energies 
so that two shocks with Mach number 10 forms later on, 
see Fig. [6] Our shock finder detects the shock fronts even 
when the actual shock is not spatially separated from the 
initial pressure jump. In the shock finder the Mach number 
is estimated via two methods: firstly it evaluates the veloc- 
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Figure 5. Shock surface area in logarithmic bins as a function 
of Mach number. Solid line depicts the measured distribution ob- 
tained from the MareNostrum simulation. Filled circles shows the 
"real" distribution, i.e., the distribution obtained by deconvolving 
the measured distribution with the data from the shock tubes. 



ity field and secondly it uses the entropy jump. In the ini- 
tial set-up the entropy jump would lead to a Mach number 
of 10, corresponding to the actual Mach number. However, 
particles needs to be accelerated to form the corresponding 
velocity field. Therefore the shock finder underestimates the 
Mach number in the beginning of the simulation. Later on, 
when the two shock fronts get close to each other, the shock 
finder includes particles from the downstream region of the 
oncoming shock, hence, it does not measure the actual up- 
stream properties anymore. A mach number estimate based 
on the velocity field would significantly overestimate the 
Mach number, since the detected velocity divergence would 
be spuriously high. However, the entropy now causes that 
the shock finder underestimates the Mach number. There- 
fore, our shock finder tends to provide only a lower estimate 
for the Mach number in complex situations, where e.g. sev- 
eral shock fronts overlap. 



6.2 Correcting the floor temperature 

One last step is to correct for the floor temperature of the 
MareNostrum simulation. Mach numbers are sensitive to 
the temperatures of the different environments. Since the 
simulation was run with a low floor temperature (200 K) 
and without any UV background, expansion cooling leads 
to very low temperatures. This can result in unrealistically 
high Mach numbers. However, the UV background heating 
leads to a simple power law between the temperature and 
the density, so the temper ature of cold regions (underdense 
regions) can be scaled by (|Hui fc Gnedinlll997h 



T 



T (l + 5) 



7-1 



(12) 
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Figure 6. Simulation of two cooling shocks. The initial densities and pressures are chosen in way that two shock fronts each with Mach 
number 10 forms when the simulation starts. 



for regions with S < 5, where 8 is the density contrast, and 
temperature To. We apply this relation to low-density par- 
ticles. Fig. [7] shows the Mach number distribution with the 
"old" floor temperature (solid lines) and with the corrected 
temperatures. Two effects can be seen: i) by correcting the 
temperature, more particles get assigned a Mach number, 
and ii) the values of the Mach numbers decreases, i.e., UV 
background heating plays an important role in cosmologi- 
ca l shocks. Nonetheless, if we compare our results to those 
of ISkillman et a l. (2008), there still high Mach numbers. In 
their simulations, they did not find shocks above M £s 200. 
This may be due to the fact that the floor temperature was 
too low, so any scaling to those low temperatures, which, 
in general, are in underdense regions, will result in a low 
increm ent of the temperature. In the case of lSkillman et all 
(2008), their floor temperature was set at T = 10 4 K, thus 
assuming the low-density gas to be ionized. 



The results in the following sections are done using the 
corrected temperature values. 



7 RESULTS 

7.1 Shock Frequency in the Cosmic Web 

We computed the surface area of identified shocks per loga- 
rithmic Mach number interval, dS(M)/ log M as described 
in Section fo.ll We divide the shock surfa ce area by t he vo l- 
ume of the simulation box, as done in iRvu et all (2003). 
The inverse of the shock surface area, 1/S, can be of 
thought as a mean separation of shocks because it is the 
simulation volume divided by the total shock surface area. 
The dS(M)/ log M distribution also indicates the frequency 
with which shocks happen in the Universe. 

Instead of studying an d classifyi n g sho cks as inter- 
nal or external as done by I Ryu et al.l (|2003T ). we use the 
Spine Web technique (see Section [3]) in order to character- 
ize voids, walls, filaments an d clusters by their morphology. 
lAragon-Calvo et al.l (|2010bT ) found that while all morpholo- 
gies occupy a roughly well-defined range in density, this is 
not sufficient to differentiate between them given their over- 
lap. Environment defined only in terms of density fails to in- 
corporate the intrinsic dynamics of each morphology. Never- 
theless, we also construct the differential shock surface area 
distribution by u sing density an d temp erature cuts in order 
to compare with lSkillman et al.l (|2008h . 

Fig.[8]shows the differential shock surface area as a func- 
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Figure 7. Histogram showing the Mach number distribution 
when the original MareNostrum floor temperature is used (solid 
line), and when applying Eqn. (dotted lines). 



tion of logarithmic Mach number bins. The top left panel 
shows the shock surface area of the different elements of the 
cosmic web as computed by the SpineWeb procedure. The 
top right panel shows the latter results after deconvolving 
them with the method outlined in Appendix [A] The bot- 
tom left panels shows the distribution divided into several 
ranges of overdensity values, while the bottom right panel 
shows the distribution divided in temperature ranges. 

We see that there is a large range of shocks of different 
strengths, with Mach numbers in the range M = 1.5 ~ 1700. 
However, it is important to notice that the large Mach num- 
bers may not occur in reality, and even the correction ap- 
plied to the temperatures in the low-density regions (see 
Section [S3] may not be enough to correctly correct the tem- 
peratures in the low density regions. 

When looking at the top panels, we observe that each el- 
ement of the cosmic web has a distinctive shock surface area 
distribution, which follows their characteristic density and 
temperature distribution and the dynamics between them. 
After deconvolving the shock surface area distribution ob- 
tained from the MareNostrum simulation using the method 
outlined in Appcndix[X] we see that there is a shift in the dis- 
tribution. In order to quantify this, we will define the char- 
acteristic Mach as the median Mach number of the different 
distributions. For example, for clusters, after deconvolving 
the shock surface area distribution with the deconvolution 
kernel, the characteristic Mach number is A4 — 1.8, slightly 
different from the value of M — 1.6 obtained without decon- 
volution. The next analysis is done using the deconvolution 
results. 

Shocks in voids have a low surface area for low Mach 
numbers in the range M m 1.5 — 10. The distribution shows 
a peak at M. » 100, extending up to M « 1300. The 
characteristic Mach number is A4 V oids ~ 18. It is expected 



that voids have strong shocks since void regions have low 
temperatures and they accrete onto walls and/or filaments, 
which are denser, hotter structures. Walls present a simi- 
lar distribution than that of voids, but they have a higher 
shock surface area for M < 10, and then decreases for larger 
Mach numbers. The characteristic Mach number of walls is 
-Mwalls » 7.5. 

On the other hand, filaments present a characteristic 
Mach number Armaments ~ 6.2, in which their surface area 
is higher than that of voids and walls. The surface area then 
decreases at M ~ 40, becoming less than that of voids and 
walls. 

Clusters have, in comparison with the rest of the cosmic 
web, weak shocks, with a maximum Mach number of M w 
10. The shocks that corresponds to the cluster environments 
are interior cluster shocks. In this region, temperatures are 
high and homogeneous, without large temperature jumps. 
The characteristic Mach number for interior cluster shocks 
is M ~ 1.8. 

We also construct the differential shock surface area in 
terms of tempe r ature and density, in order to compare with 
Skill man et all (|2008T ). This is shown in the lower panels of 
Fig. [8] The difference between these panels and the upper 
panels is that here we divide the shock surface area in tem- 
perature and density ranges. The temperature and density 
ranges are intended to denote the different temperature and 
density of the elements of the cosmic web. for simplicity, 
we do not deconvolve these distributions. When looking at 
density cuts (lower left panel of Fig. [HJ, it is interesting to 
notice that cluster cores (5 > 10 4 ) have low Mach numbers 
and low surface area. This is mainly due to the fact that 
this region has already been heated, and so its temperature 
is high. Low-density regions (5 < 1) may be associated with 
voids. These regions have high characteristic Mach numbers 
due to their low density and low temperature. We see that 
there is an increase of the characteristic Mach number as 
the density decreases. 

Similar results are observed when looking at the tem- 
perature cuts. Each cut has a characteristic Mach num- 
ber, which increases as temperature decreases. For high 
temperature (T > 10 7 K), the characteristic Mach num- 
ber is A4 char » 2, w hile for T < 1000 K, AW « 20. 
ISkillman et al.l (|2008h states that this characteristic Mach 
number is due to the maximum temperature jump possible 
for a given temperature. 

In general , we fin d good agreement with the result of 
ISkillman et al.l ()2008h for cuts with low density and tem- 
perature. For dense structures we find tha t our results dif- 
fer significantly from ISkillman et alj |2008l ). More precisely, 
we find that even in regions with 6 > 10 3 s hocks with 
M > 100 can be found, while ISkillman" et al l l|2008l ) find 
that shocks have Mach numbers basically below 10 for re- 
gions with 8 > 10 3 . Comparing the density cuts (lower left 
panel of Fig. [8} with the upper panels we see that dense re- 
gions with high Mach number does not reside in clusters. In- 
stead, these are dense structures in filaments, walls or voids. 
The particle based TreeSPH simulation technique, as used 
in Gadget, tends to from clumpy structures. Therefore, we 
find dense structures with 5 > 10 3 even outside of clusters. If 
these structures reside close to an accretion shock, we find 
a high Mach number at a relatively high density. This ef- 
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Figure 8. Differential shock surface area as a function of logarithmic Mach number bins. The top left panel shows the distribution for the 
different elements of the cosmic web: voids (dash-dotted line), walls (dashed line), filaments (dotted line) and clusters (dash-dot-dotted 
line). The solid line shows the combine distribution of all the elements of the cosmic web. By applying the deconvolution technique 
(Sect. [Sj) to these distributions, we obtained the results presented in the top right panel. The bottom left panel shows the shock surface 
distribution divided into several ranges of overdensity. The bottom right panel shows the distribution divided into several temperature 
ranges. 



feet is enhanced by smearing out shock fronts in SPH by 
approximately two times the smoothing length. 

For reg ion with 5 > 10 4 w e find significantly less 
shocks than ISkillman et alJ |2008). The difference is likely 
cause d by a lower e f fective resolution in our SPH simula- 
tion. ISkillman et al.l (|2008T ) used a simulation box with a 
side length of 512 h^ 1 Mpc. They used a root grid resolu- 
tion for both dark mater and gas of 512. For dark matter 
the MareNostrum simulation is significantly better resolved. 
In the AMR scheme with a density refinement criterium 
the resolution follows the mass roughly in a similar way as 
SPH. However, since the SPH kernel smoothes over N par- 
ticles (64 in our case), the effective resolution for gas is in 
the MareNostrum simulation by a factor of 2 low e r than in 
the 'Santa Fe Light Cone' used by ISkillman et all l |200Sf ). A 
higher resolution as used in the MareNostrum simulation is 
necessary to form shock fronts in regions with 6 > 10 4 . 



7.2 Temperature Distribution and Physical 
Properties of the Cosmic Web 

As mentioned earlier, one of the big advantages of the 
SpineWeb procedure is that we can precisely characterize 
the different elements of the cosmic web. We can, therefore, 
study the temperature and density distribution in the dif- 
ferent elements of the cosmic web. 

Table [2] shows the mean values of the temperature and 
density contrast of the MareNostrum simulation. The shown 
temperature is the simulation temperature corrected using 
the relation in Eqn. [T2j The mean values were obtained by 
averaging the temperature and density of the particles that 
belong to the different environments. As expected, clusters 
are the hottest structures in the Universe, with a mean tem- 
perature of 7.90 x 10 7 K. This is because clusters have deep 
potential wells into which baryons accrete, thus heating the 
clusters. 

The most salient features of the cosmic web are the large 
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Table 2. Mean temperature and mean density contrast of the 
different elements of the cosmic web at z = 0. 



<T> (K) a T (8) a s 



Clusters 7.90xl0 7 8.13xl0 7 885.73 1108.28 

Filaments 7.30xl0 6 1.44xl0 7 71.87 180.27 

Walls 2.21X10 6 6.84X10 6 18.71 71.24 

Voids 6.96x10 s 4.30xl0 6 2.89 25.17 
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Figure 9. Temperature distribution of the different environments 
of the Cosmic Web. 



filamentary networks which are interconnected across vast 
distances. The filamentary network permeates all regions of 
space, even the underdense voids. We find that filaments 
have a mean temperature of 7.30 x 10 6 K, which is well within 
the range estimated for the filamentary WHIM temperature 
(10 5 - 10 7 K). 

Table [5] also shows the density contrast of the con- 
stituents of the cosmic web. Filaments are the second densest 
objects in the Universe, with (8) ~ 70. Walls are less dense 
than filaments, with (8) ~ 19. Finally, voids have a mean 
density contrast of (5) ~ 3. 

Fig. [13] depicts a zoom-in of box number 3 of Fig. 1101 
As in the previous figure, we only consider void parti- 
cles as defined by the Spine Web procedure. Interesting is 
to notice that voids a re not empty. As mentioned before, 
iGottlober et~ai1 (|2003l ) stated that voids are like miniature 
universes. Voids have what may be called sub- walls and sub- 
filaments (a sub-network), with void galaxies the nodes of 
the sub-filaments. These sub-walls and sub-filaments do not 
have a high enough density to be considered walls or fila- 
ments on their own. In the region chosen, of 60h~ Mpc on 
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Figure 10. Density field of the MareNostrum Universe. The fig- 
ures represents a slice of 20/i~ 1 Mpc thick through the z-axis. 
Three elements are distinguished: the most massive cluster in the 
simulation (box number 1), a region of rich filaments (box number 
2) and a void region (box number 3). The different environments 
were find using the SpineWeb technique. 
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Figure 11. Temperature and Mach number in an infinitcsimally 
slice for the most massive cluster in the MareNostrum simulation 
(box 1 in Fig. [T0t. 
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Figure 12. Temperature and Mach number in a thin slice of a 
rich region of filaments in the MareNostrum simulation (box 2 in 
Fig. [Toll 



a side, we see what it could be void galaxies and some tenu- 
ous filamentary structure. The temperature of this region is 
T » x 10 3 — 10 6 . We see that a void region has high Mach 
numbers mainly due to the low temperature of the accreting 
gas. We find that void regions in the MareNostrum simula- 
tion have shocks with Mach numbers M > 1000. However, 
these high Mach numbers may not be real, and they may be 
due to the low temperature floor set in the simulation (see 
Section [7UJ. 

7.3 Radio Emission in the Cosmic Web 

Each SPH particle has been assigned a radio luminosity via 
Eqn. [5] and with the correct setting of /a(-M) (Eqn. [TTJ. 
With this in hand it is possible to construct artificial radio 
maps. To this end, we project the emission of each parti- 
cle long the line of sight and smooth it with the SPH ker- 
nel size. The radio emission is computed for an observing 
frequency of 1.4 GHz and a hypothetical beam of 10x10 
arcsec 2 . For comparison, we also compute contours of the 
bolometri c X-ray flux. Each S PH particle is assigned a lu- 
minosity (N avarro et al.lll995T ) 
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Figure 13. Temperature and Mach number in a 10h~ x Mpc thick 
slice for a void region (box 3 in Fig. I10p. 

L X = 1.2X10 ergs — ) . (13) 

We will study the radio emission in each component of 
the cosmic web. The results in these section are obtained 
using an electron efficiency of £ c = 0.005. In Section \7A\ we 
calibrate our simulated data against observations, yielding a 
different electron efficiency and, therefore, different number 
of radio objects. 

7.3.1 Clusters 

Of the total sample of galaxy clusters present in the 
MareNostrum Universe (~ 3800) we find that only 565 clus- 
ter of galaxies have radio objects (693) with Pi. 4 > 10 30 ergs 
s _1 Hz -1 (see Sec. for a correction on this number based on 
a correction of £ c ). This represents only ~ 15% of the total 
galaxy cluster sample. In total, we find that ~ 41% of the 
total sample of galaxy clusters host diffuse radio emission. 
These clusters are not necessarily the most massive ones, 
but they are within the entire mass range. 

Fig. [14] shows the synthetic observation of X-ray and 
radio emission of the most massive cluster in the simulation. 
The thin straight lines in the maps indicate the direction of 
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Figure 14. Synthetic observations of X-ray and radio emission. 
The radio flux at 1.4 Ghz is computed for a hypothetical beam 
size of 10 X 10 arcsec 2 . Upper panel: the galaxy cluster depicted 
in Fig. 1111 Middle panel: the filamentary region shown in Fig. 1121 
Bottom panel: the void region shown in Fig. 1131 The thin straight 
lines indicate the direction of polarization. The length of the line 
indicates the degree of polarization. 




Cluster X-ray temperature [keV] 

Figure 15. Radio luminosity of radio objects versus the emission 
weighted temperature of galaxy clusters. 

polarization. The length of the line indicates the degree of 
polarization. For comparison the length corresponding to 
10%, 50% and 100% is shown. The polariza tion is comp uted 
according to the f ormalism described in (|Burnlll96d ) and 
ijEnfilin et alJll998h . The radio emission is caused by internal 
shock fronts in the clusters. We could, therefore, exclude 
external shocks as sources for prod ucing radio emission in 
cluster (see also iHoeft et alj l200ct ) . We find that clusters 
host the majority of radio objects, with ~ 20 — 30 objects 
with P1.4 > 10 32 ergs s" 1 Hz -1 (see Section 171} . Of a total 
of ~ 1500 radio objects with L R > P XA > 10 30 ergs s" 1 
Hz" 1 found in the simulation, ~ 500 are located in galaxy 
clusters, i.e., 33% of the total radio objects. 

We also estimate the relation between the radio lumi- 
nosity, Pi. 4, of the radio objects and the emission weighted 
temperature of its galaxy clusters hosts. This r esult should 
be similar to the one obtained bv lHoeft et al. I (|2008h . only 
that our sample of galaxy clusters is 10 times larger than 
theirs. 

Fig. \T5\ shows that those clusters that host the very 
luminous radio objects, i.e., Pi. 4 > 10 32 ergs s _1 Hz -1 , are 
very hot (Tx > 8 keV)., i.e . , very massive. As was also 
pointed out by iFeretti et al.l l|2004l ). there is a correlation 
between the radio luminosity of the radio objects and the 
temperature of its host cluster, i.e., the hottest clusters hosts 
the most luminous radio objects. 

Our re s ults a re quite similar to those presented in 
IHoeft et al.l l|2008l ). This is expected since their sample is 
a fraction of the sample of clusters presented here. The re- 
sults presented here reinforce that our simple emission model 
recreate diffuse radio emission not only in galaxy clusters, 
but also in the cosmic web, and picks up the trend in the 
radio luminosity - X-ray relation. 

7.3.2 Filaments 

iBagchi et~a l. (2002) reported the existence of a large-scale 
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diffus e radio em i ssion from a large-scale filamentary net- 
work. iKim et al.l (|l989T ) also suggested the presence of dif- 
fuse radio emission in the Coma Supercluster. Temperatures 
of filaments are in the range T « 10 5 - 10 7 , and magne- 
tohydrodynamic simulations (MHD) esti mate a magnetic 
field that varies fro m 10" 1 nG to 10 n G llSigl et al 120031 : 
iBriiggen et"alll2005f ). This is in line with iRvu et all <|2008h . 
who suggest an aver age magnetic streng th to be of the order 
of 10 nG. However, iBaechi et aD (|2002h estimated that the 
strength of the magnetic field in the filamentary network 
isBfs 0.3 - 0.5 j«G. In the same line, iDolag et all (|2005h 
suggests the strength of the magnetic field of filaments that 
connect galaxy clusters to be as large as 0.1 fiG. The middle 
panel of Fig. [14] shows the synthetic observations of the X- 
ray and the radio emission for the filamentary region shown 
in Fig. 1121 We found that filaments show diffuse radio emis- 
sion, especially in areas near to cluster of galaxies. However, 
although filaments have strong Mach numbers (see Fig. I12|) . 
the flux in "pure" filamentary regions, i.e., regions far from 
the outskirts of clusters is low, as seen in Fig. 1141 and may 
not be detectable with the upcoming radio telescopes. 

As for the number of radio objects, we find that fila- 
ments host fewer radio objects than clusters, with ~ 2 — 3 ra- 
dio objects with P 1A > 10 32 ergs s _1 Hz -1 (see Section PO)) . 



7.3.3 Voids 

As discussed in Sect. [5] esitmation of magnetic field out- 
side galaxy clusters is difficult. Moreoever, MHD simulations 
have not been conclusive as well, assi gning a magneti c field 
that varies from 10~ 3 nG to lO^ 1 nG (|Sigl et aljr2003h . Fur- 
ther studies on the magnetic field on voids are necessary to 
precisely determine the strength of these magnetic fields. 

As can be seen in the bottom panel of Fig. 1141 there 
are some areas in the void region that shows a very tenuous 
radio emission, with fluxes < 10 -4 mjy. This is basically due 
to two facts: temperatures in voids, as we have seen, is very 
low compared to the other components of the cosmic web. 
It has also been suggested that the magnetic field in voids 
is quite low. Although observationally it has been possibly 
to estimate the mag netic field in clusters of galaxies (e.g., 
iBonafede et al.l 2010), this has not been possible for voids. 



7.4 Radio objects luminosity function 

As we have shown, there exists the possibility of finding 
radio emission in other environments of the cosmic web, such 
as filaments. We now turn our attention to the possibility 
of finding individual radio objec ts in the cosmic we b. We do 
this in a similar way to that of iHoeft et al.l (|2008l '). with a 
slight difference. We select all particles whose radio emission 
lies above a very low threshold. This threshold corresponds 
roughly to particles with M > 2. Instead of finding their 
nearest neighbor a nd link them accord ing to their smoothing 
length, as done bv lHoeft et ail <|200Sl) . we use HOP to link 
the particles and find individual groups. This is a somewhat 
simpler approach, but given that HOP works on a density 
basis, it assures us that the particles will form a compact 
group. Having found the groups, we compute the cumulative 
number of radio objects above a given luminosity, that is, 
the luminosity function of diffuse radio objects. 
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Figure 16. Cumulative number density of radio objects. The 
thick solid line indicates the total number of radio objects in the 
simulation box. The thin solid line shows the number of radio ob- 
jects within clusters. The dashed line indicates the radio objects 
present in filaments. The dashed-dotted line depicts the radio ob- 
jects in walls, while the dashed-dotted-dashed line represent the 
radio objects present in voids. For comparison, the NVSS relic 
luminosity function (dotted line). 



Fig. [16] shows the radio luminosity function for radio 
objects in the simulation box. For comparison, we also show 
a radio relic luminosity function constructed with relic 



data observations (Giovannini et al.l 19991; Govoni et al.l 


|200ll: Bonafede et all 120091: 


Clarke & Ensslinl 120061: 


iFeretti et al.l 120061: iBagchi et al. 


20061: iRottgering et all 


19971 van Weeren et al.l 12009 


: Feretti et all 200ll: 


Giacintucci et all I200SL IVenturi et all l2007t). Most of 



these dat a were obtained fro m the NRAO VLA Sky Survey 
(NVSS) i|Condon et alJll998h . 

We see that with our simple assumptions, we find a 
much larger number of radio objects in the entire simula- 
tion box. Also, if we compare with our cluster radio object 
luminosity function (thin solid line), we still find more radio 
objects. Galaxy clusters host the majority and the most lu- 
minous objects, with ~ 34 objects with Pi. 4 > 10 32 ergs s _1 
Hz -1 . Filaments also host some radio objects, but most of 
them are less luminous, with only a couple with Pi. 4 > 10 32 
ergs s _1 Hz -1 . Finally, walls and voids also host a few radio 
objects, but they are rather weak. Table [3] shows the frac- 
tion of radio objects by environment for different total radio 
luminosity ranges. It can be seen that clusters host the ma- 
jority of more luminous radio objects (Pi. 4 > 10 32 ergs s _1 
Hz -1 ), but when considering the complete sample of radio 
objects, these are located within filaments. This is to expect, 
since a great part of the mass of th e Universe resides in fil- 



aments ( Aragon-Calvo et al.ll2010bh . This is also an indica- 



tion that the detection of diffuse radio emission in filaments 
could shed light on the missing-baryon problem. Walls also 
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Figure 17. Same as Fig, 1161 but the total radio power has been 
rescaled by using a electron efficiency parameter of e e = 0.0005. 
By doing this, radio objects are less luminous, and we are able to 
match our simulated relics with observational NVSS relics. 



host a great number of low-luminous radio objects, while in 
voids the presence of r adio objects i s neg ligible. 

As pointed out bv lHoeft et al.l l|200ct ). the overestima- 
tion of radio objects may be due to an overestimation of 
the electron efficiency (£ c ), as well as the magnetic fields. 
The value of the the electron efficiency parameter is quite 
uncerta in from current obse rvational and theoretical con- 
strains (|Skillman et alJfeoiOr ). From Eqn. [8] we see that the 
relationship between the total radio power and £ c is linear, 
so in order to obtain less luminous objects, we just have to 
rescaled the total radio power, e.g., an electron efficiency of 
£o = 0.0005 will lead to a total radio power that is low by a 
factor of fO. In order to match our results to observations, 
i.e., the NVSS radio relic luminosity function, we rescale 
the luminosity function using e e — 0.0005, shown in Fig. 1171 
As expected, radio objects are less luminous, and therefore, 
less radio objects are found. In the plotted range, there are 
no radio objects in voids and filaments, only 1 in filaments, 
and some group of ten in clusters. With this electron effi- 
ciency, the cluster radio objects matches the observational 
NVSS radio relic luminosity function (see Fig. ll7[l . However, 
we predict more relics at the high and low power compared 
with the observed function. This can be to due to the fact 
that the observed relics, to date, is still small, compared 
to the large number of relics found in the simulation with 
our method. Upcoming radio observations will increase the 
number of observed relics, making the statistics stronger. 
Table [3] shows the fraction of radio objects by environment 
for different total radio luminosity ranges with the adjusted 
electron efficiency parameter. A study of the effects of the 
magnetic fields is beyond the scope of this paper, so we do 
not investigate this. 

The luminosity function (Eqn. 1 16 p suggest what can be 



Table 3. Fraction of total radio objects (~ 14600) with a given 
radio luminosity Pi. 4 in terms of their environment assuming an 
electron efficiency of £ c = 0.005 (£ c = 0.0005) Units of radio 
luminosity given in ergs s -1 Hz . 



> 10 3 



10 30 < P1.4 < 10 3 



Pi. 4 < 10 3 



Clusters 
Filaments 
Walls 
Voids 



0.23% (0.03%) 
0.01% (0%) 
0% (0%) 
0% (0%) 



4.51% (1.6%) 
3.65% (0.3%) 
0.18% (0%) 
0.03% (0%) 



10.33% (13.46%) 
74.24% (77.58%) 
6.40% (6.58%) 
0.42% (0.45%) 



expected with the use of the upcoming radio telescopes. An 
increase of the surface brightness sensitivity by a factor of 10 
will result in an increase of the number of radi o objects by a 
factor of 10-100 (see also lSkillman et al.l l|2010t ). We can also 
estimate the number of radio objects to be found in a given 
survey area and redshift depth by solving the equation: 



dV 
dzdQ 



(*) = D 1 



(l + z) 2 D 2 , 
E(z) 



(14) 



where Dh = c/Ho is the Hubble distance, D a 
is the angular diameter distance and E(z) = 
\/Q.m{l + z) 3 + + z) 2 + Qa. For a flat, A-dominated 
Universe like the MareNost rum simulation, fl k = and 
Eqn. [2] can be rewritten as (|Carroll et al.lll992p 
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(15) 



where Dm is the comoving distance, which, for a flat Uni- 
verse, coincides with the line-of-sight comoving distance. 
Eqn. [^corresponds to the total comoving volume in an all- 
sky survey out to redshift z. For the MareNostrum Universe, 
the result of Eqn. [14] corresponds to an all-sky survey out to 
z ~ 0.18. If we assume an electron efficiency of £ c = 0.005 
(£ e = 0.0005), we expect to find 34 (4) radio objects in 
30 (3) galaxy clusters, 2 (1) radio objects in filaments and 
none in walls and voids, all with a total radio luminosity 
of Pi. 4 > 10 32 ergs s _1 Hz -1 within this cosmological vol- 
ume. Assuming that the number den sity of radio objects is 
nearly constant through cosmic time (|Skillman et al.ll2010l ). 
we could estimate that in an all-sky survey out to z — 0.5 
(which corresponds to a volume of ~ 15 times the previous 
one), we could find 510 (60) radio objects in 450 (45) clus- 
ters, 30 (15) radio objects in filaments and a few in walls 
and voids. 

The use of the Spine Web enable us to make predictions 
of the radio flux in filaments. Assuming an spectral index of 
a = 1, filaments at a redshift z ~ 0.15, in a frequency of 150 



MHz, a radio flux of S15 
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0.12 pjy. 



We have analyzed one of the largest hydrodynamical simu- 
lations of cosmic structure formation, namely the MareNos- 
trum simulation, with the purpose of studying the possibil- 
ity of detecting radio emission in the cosmic web. To that 
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purpose, we use: i) a novel method for estimating the ra- 
dio emission of strong shocks that occur during the pro- 
cess of structure f ormat ion in the Universe developed by 
iHoeft fc Briiggenl |2007n and ii) the SpineWeb technique 
The SpineWeb procedure deals with the topology of the un- 
derlying density field, correctly identifying voids, walls and 
filaments. While previous studies focused on the study of 
cosmological shocks in selected catalogs of galaxy clusters 
or by defining the different environments of the cosmic web 
by density and/or temperature cuts, we are able to correctly 
disentangle the complex filamentary network of the cosmic 
web. The computation of the radio emission is based on an 
estimate for the shock surface area. The radio emission is 
computed per surface area element using the Mach number 
of the shock and the downstream plasma properties. 

In order to properly calculate the radio emission, two 
corrections were made. The first correction deals with the 
shock surface area. The constant area factor (f a ) of Eqn. [8] 
was corrected using shock tube tests and a convolution 
method described in Appendix [X] One deals with the tem- 
perature of the simulation. Given that the simulation was 
run with a low floor temperature and without UV back- 
ground, expansion cooling leads to very low temperature, 
which leads to very high Mach numbers. We corrected this 
by scaling the temper ature of cold, low den sity regions us- 
ing T = T (l + 5) 7 " 1 l|Hui fc Gnedinlll997r ). where 5 is the 
overdensity. The first step was to calculate the strength 
of cosmological shocks, characteriz ed by their Mach num- 
ber using the method described in IHoeft et al.l (|2008l ). We 
find that each environment of the cosmic web has a dis- 
tinct Mach distribution: voids have a characteristic Mach 
number of Avoids ~ 18; walls, A^waiis ~ 7.5; filaments, 
A^Hiamcnts ~ 6.2; and clusters, A^dusters ~ 1-8. The strength 
of the shock, i.e., the characteristic Mach number, is closely 
related to the temperature and density of the medium. Voids 
are the coldest and less dense region in the Universe, with a 
mean temperature (T) vo id s ~ 7 x 10 s K and (5} V oids ~ 2.9, 
and they interact with walls and filaments, which are denser 
and hotter structures, with filaments reaching a mean tem- 
perature of (T) voids « 7.3 x 10 6 K and (<5) vo id s « 72. 

Applying a radio emission model based on DSA to the 
simulation we are able to reproduce the diffuse radio emis- 
sion in gala xy clusters using a lower electron efficiency to 
that used bv lHoeft et all (|2008l ). This diffuse radio emission 
is also present in filaments. In walls and voids the emission 
is very weak, resulting in almost no presence of radio ob- 
jects. By identifying radio objects based on their density, 
we find that clusters hosts the majority of such objects. As- 
suming an electron efficiency of £ c = 0.005, we find that 
30 galaxy clusters (out of 3865) host 34 radio objects with 
Pi. 4 > 10 32 ergs s _1 Hz -1 . If we consider a radio luminos- 
ity of 10 30 < Pi. 4 «S 10 32 ergs s" 1 Hz" 1 , then 659 radio 
objects are found in 547 galaxy clusters, only ~ 15% of 
the total sample. These radio objects are present in the en- 
tire mass range of the galaxy clusters, i.e., both low mass 
and high mass clusters present radio objects. Filaments also 
host radio objects, although only 2 of them are very lumi- 
nous (Pi. 4 > 10 32 ergs s _1 Hz -1 ). However, the majority 
of low luminous radio objects (Pi. 4 < 10 30 ) are present in 
filaments. The presence of radio objects in walls and voids 
is rather negligible. 

We calibrated our simulated data against the observa- 



tional NVSS radio relic luminosity function, yielding an elec- 
tron efficiency of £ e = 0.0005. Th i s is a factor of 10 less than 
the one estimated bv lHoeft et al.l (|2008l ). When using this ef- 
ficiency, the number of very luminous radio objects drop dra- 
matically. We only find 4 radio object with Pi. 4 > 10 32 ergs 
s" 1 Hz" 1 and 1 in filaments. In the range 10 30 < P1.4 sC 10 32 
ergs s _1 Hz -1 , we find 235 radio objects in 198 galaxy clus- 
ters (5% of the total sample) and 43 in filaments. 

An all-sky survey up to z — 0.5 and assuming an ef- 
ficiency of £ c = 0.005 ( £c = 0.0005) should result in the 
discovery of 510 (60) radio objects in 450 (45) clusters, and 
30 (15) radio objects in filaments. With the increase of sur- 
face brightness sensitivity of the upcoming radio telescopes, 
the detection of radio objects should increase by a factor of 
10, which opens in the possibility of finding a considerable 
number of radio objects in filaments. Furthermore, we pre- 
dict that the radio flux of filaments at redshift 2 ~ 0.15, and 
at a frequency of 150 MHz, should be SibomBz ~ 0.12 /ijy. 
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APPENDIX A: CORRECTING THE SHOCK 
SURFACE AREA DISTRIBUTION 

Convolution is a mathematical operation on two functions, 
/ and g, which produces a third function that is a modified 
version of the one of the original functions. Usually, one of 
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the two functions is taken to be a kernel function which acts 
on the other function, modifying it. The standard convolu- 
tion states that 

/oo 
f(r)g(t-r)dr, (Al) 
- oo 

where * is the convolution operand. The inverse of this op- 
eration is called deconvolution. In our case, we have a mea- 
sured distribution, namely dS / d\ogM, obtained from the 
simulation and we wish to obtain the real distribution, i.e., 
the shock surface area distribution taking into account the 
particles that were assign a low Mach number in strong 
shocks and the underestimation of the different Mach num- 
bers. Therefore, we need to find the shock tube distribu- 
tion kernel taken from the shock tubes data in order to de- 
convolve it with our measure function. However, instead of 
finding the inverse of the kernel function, we will iterate the 
convolution procedure, which will yield the real distribution. 
The convolution we want to solve is 

r oo 

q m (logM) = / dlogM'f(logM,logM')q r {logM'),(A2) 
Jo 

where q rn (logM) is the measured shock surface distribu- 
tion, /(log M, log M') is the kernel function and g r (log.M) 
is the real distribution we want to find. In Eqn. IA2l we have 
taken two things into account: i) there are not negative Mach 
numbers, and ii) the kernel function depends on the Mach 
number. 

In Section [6] we saw that a simple representation to 
the shock surface area of the shock tubes are a Heaviside 
function and a Gaussian. Therefore, our kernel has the form 

Y{\ogM) lJ = K ■ H{\ogMi) + (log .Mi, log M 3 ) , (A3) 

where Qij is a Gaussian with form 

Q.. _ e -(logA1-logAl') 2 /2<T(logAl') 2 (A4) 

where the standard deviation also depends on the measured 
Mach numbers, i.e., a = cr(log.M) as seen in the previous 
section. 

After constructing the kernel with the fitting parame- 
ters of the various shock tubes, we give a guess function in 
order to start the iteration procedure. As a first guess, we 
choose q[ = q m , and convolve that with the kernel, giving as 
a result . The new guess is then calculated by subtracting 
the previous guess with the difference between the result of 
the convolution of the first guess with the measured func- 
tion, i.e., q' 2 — q\ — (g™ — q m )- We found that after iterating 
5 times we are able to recover q m . 
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